
      !<@brief Gauss Newton Method
      subroutine GN(N,X,tol,itmax,Hnnz)

      implicit none
      integer :: N !< [in] número de variáveis de otimização
      double precision :: X(N) !< [in,out] vetor de parâmetros de otimização

      double precision:: tol !< [in] critério de parada, tolerância para o valor da função objetivo
      double precision:: itmax !< [in] critério de parada, número máximo de iterações

      double precision, dimension(:), allocatable:: hdd !< direção de descida

      double precision :: gradF(N) !<  gradiente da função objetivo em relação às variáveis de otimização
      integer :: Hnnz !<  número de elementos não nulos da Hessiana
      double precision :: Hvalue(Hnnz)!<  Elementos não nulos da aproximação da Hessiana para ser usada como matriz de precondicionamento
      integer :: Hrowind(Hnnz) !<  Índice das linhas dos elementos não nulos da Hessiana aproximada. Esquema CCS de armazenamento.    
      integer :: Hcolptr(N+1) !<  Guarda o nnz do primeiro elemento não nulo de cada coluna. Esquema CCS de armazenamento. 

      double precision, dimension(:), allocatable :: Xnew !< vetor de parâmetros de otimização atualizado
      double precision:: gamma
      double precision:: VALUE
      double precision :: fobj
      double precision :: fobjnew
      integer :: it=0

      integer:: i,k
      integer:: aval



      integer :: iopt 
      integer :: nrhs
      integer :: ldh
      integer :: info
      integer*8 :: factors



      allocate(hdd(N))

      allocate(Xnew(N))

      fobj = 100.0d0

      do
	  if(it.gt.itmax .or. fobj.lt.tol)exit

	  call BOTHGN(fobj,gradF,X,it)

	  do i=1,N
	    hdd(i) = - gradF(i)
	  enddo

	  ldh = N
	  nrhs = 1

	  iopt = 1! fatorar
	  call c_fortran_dgssv( iopt, N, Hnnz, nrhs, Hvalue, Hrowind, Hcolptr, hdd, ldh, factors, info) !CCS  --- no c_fortran_zgssv estou fazendo H^T . Y = Z, mas nao tem problema porque H eh simetrica
	  if (info .eq. 0) then
	      write (*,*) 'GN - H Solve succeeded'
	  else
	      write(*,*) 'GN - H INFO from factorization = ', info
	      stop
	  endif

	  iopt = 2! resolver
	  call c_fortran_dgssv( iopt, N, Hnnz, nrhs, Hvalue, Hrowind, Hcolptr, hdd, ldh, factors, info) !CCS  --- no c_fortran_zgssv estou fazendo H^T . Y = Z, mas nao tem problema porque H eh simetrica
	  if (info .eq. 0) then
	      write (*,*) 'GN - H Solve succeeded'
	  else
	      write(*,*) 'GN - H INFO from factorization = ', info
	      stop
	  endif

	  iopt = 3! limpar
	  call c_fortran_dgssv( iopt, N, Hnnz, nrhs, Hvalue, Hrowind, Hcolptr, hdd, ldh, factors, info) !CCS  --- no c_fortran_zgssv estou fazendo H^T . Y = Z, mas nao tem problema porque H eh simetrica
	  if (info .eq. 0) then
	      write (*,*) 'GN - H Solve succeeded'
	  else
	      write(*,*) 'GN - H INFO from factorization = ', info
	      stop
	  endif
 	      write(*,*) 'abntes*****',aval,gamma, fobj, fobjnew 
	  gamma = 1.0d0
	  aval=0
	  do 
	      do i=1,N
		  Xnew(i) = X(i) + gamma*hdd(i)
	      enddo	
	      fobjnew=VALUE(Xnew)
	      write(*,*) '*****',aval,gamma, fobj, fobjnew
	      if(fobjnew.lt.fobj .or. aval.gt.10)exit
	      aval= aval+1
	  enddo

	  do i=1,N
	      X(i) = Xnew(i)
	  enddo	

	  it = it + 1
      enddo


      deallocate(hdd)
      deallocate(Xnew)

      return
      end subroutine GN